/*!  \file  prototypes.h
 *	  \brief  contains all blas lapack prototypes used in this project
 * 
 *   \details   details
 *   \author    Yordan Madzhunkov German Research School for Simulation Sciences GmbH
 *   \version   1.0
 *   \date      03/29/2013 05:47:13 PM
 *   \copyright GNU Public License.
 */
#pragma once



/************
 *  B L A S *
 ************/

extern double ddot_( const int *N, const double *a, const int *inca, const double *b, const int *incb );
extern void  dgemm_(char*, char*, int*, int*, int*, double*, double*, int*, double*, int*, double*, double*, int*);

/***************
 * L A P A C K *
 ***************/

// computes the svd factorization of general rectangular matrix
//extern void dgesvd_( char* jobu, char* jobvt, int* m, int* n, double* a,
//			int* lda, double* s, double* u, int* ldu, double* vt,int* ldvt,
//			double* work, int* lwork, int* info );


//  DSTEV computes all eigenvalues/eigenvectors of a   
//      real symmetric tridiagonal matrix A. 
//extern int dstev_(char *jobz, int * n, double * d, double *e, double * z, 
//	int* ldz, double  * work, int *info) 
// DSTEVD computes all eigenvalues/vectors of a real symmetric tridiagonal matrix. Uses a divide and conquer algorithm.

// SVD
extern void dgesvd_ ( //computes SVD MxN matrix A, A = U * SIGMA * transpose(V)
char   * jobu , // 'A': All columns, 'S'  the first min(m,n) of U/V are returned in array U/V
char   * jobvt, // 'O':  the first min(m,n) rows of U/V^T, 'N' : no singular vectors
int    *     m, // Rows of the input matrix A
int    *     n, // Columns of the input matrix A
double *     A, // (input/output) On entry, the M-by-N matrix A., at output content is destroyed
int    *   ldA, // leading dimension of the array A, should be M
double *     s, // double[min(M,N)] containg singular values
double *     u, // left singular vectors, stored columnwise
int    *   ldU, // 
double *    vt, // right singular vectors, stored rowwise
int    *  ldvt, // 
double *  work, // working array of size lwork
int    * lwork, // The dimension of the array WORK. LWORK = -1,then a workspace query is assumed
int    *  info  //  0  :  successful exit.
                // < 0 :  if INFO = -i, the i-th argument had an illegal value.
                // > 0:  if DBDSQR did not converge, INFO specifies how many
                // superdiagonals of an intermediate bidiagonal form B
                // did not converge to zero. See the description of WORK
                // above for details. 
);
extern void dgejsv_(// preconditioned Jacobi SVD algorithm. 
char * joba, // 'C' : This option works well (high relative accuracy)
 				 // 'E' : Computation as with 'C' with an additional estimate of the realistic error bound.
       		 // 'F' : If input matrix structure is not known, and relative accuracy is desirable, then this option is advisable.
       	    // 'G' : Computation as with 'F' with an additional estimate of the too pessimistic error bound.
       		 // 'A' : Small singular values are the noise and the matrix is treated as numerically rank defficient. The error in the computed
				 // singular values is bounded by f(m,n)*epsilon*||A||. The computed SVD A = U * S * V^t restores A up to f(m,n)*epsilon*||A||.
       		 // 'R': Similar as in 'iA'. Rank revealing property of the initial  QR factorization is used do reveal (using triangular factor)
             // a gap sigma_{r+1} < epsilon * sigma_r in which case the numerical RANK is declared to be r. The SVD is computed with absolute 
             // error bounds, but more accurately than with 'A'.
char   * jobu , // 'U': N columns of U are returned in the array U.
       			 // 'F': full set of M left sing. vectors is returned in the array U.
       			 // 'W': U may be used as workspace of length M*N. See the description of U.
       			 // 'N': U is not computed.

char   * jobvt, // 'V': N columns of V are returned in the array V; Jacobi rotations are not explicitly accumulated.
       			 // 'J': N columns of V are returned in the array V, but they are computed as the product of Jacobi rotations. Allowed only if JOBU != 'N'
					 // 'W': V may be used as workspace of length N*N. See the description of V.
       			 // 'N': V is not computed.

char   * jobz , // 'A': All columns, 'S'  the first min(m,n) of U/V are returned in array U/V
					 // 'O':  the first min(m,n) rows of U/V^T, 'N' : no singular vectors
int    *     m, // Rows of the input matrix A
int    *     n, // Columns of the input matrix A
double *     A, // (input/output) On entry, the M-by-N matrix A., at output content is destroyed
int    *   ldA, // leading dimension of the array A, should be M
double *     s, // double[min(M,N)] containg singular values
double *     u, // left singular vectors, stored columnwise
int    *   ldU, // 
double *    vt, // right singular vectors, stored rowwise
int    *  ldvt, // 
double *  work, // working array of size lwork
int    * lwork, // The dimension of the array WORK. LWORK = -1,then a workspace query is assumed
int    * iwork, // INTEGER array, dimension  dimension M+3*N
int    *  info  //  0  :  successful exit.
);
/*  TO DO

extern void dggsvd_ ( 
char * jobu, 
char * jobv, 
char * jobq, 
int * m, 
int * n, p, k, l,
double * A, 
int * ldA, 
B, 
int * ldB, 
alpha, 
beta, 
u, int * ldU, v, int * ldV,
Q, int * ldQ, double * work, iwork, int * info )

SVD, divide-and-conquer
dgesdd_ ( char * jobz, int * m, int * n, double * A, int * ldA, s, u, int * ldU, vt, int * ldvt, double * work, int * work, iwork, int * info )

SVD, Jacobi
dgejsv_ ( char * joba, char * jobu, char * jobv, char * jobr, char * jobt, char * jobp, int * m, int * n, double * A, int * ldA, sva, u, int * ldU, v, int * ldV, double * work, int * work, iwork, int * info )
dgesvj_ ( char * joba, char * jobu, char * jobv, int * m, int * n, double * A, int * ldA, sva, mv, v, int * ldV, double * work, int * work, int * info )


Reduce to bi-diagonal
dgebrd_ (       int * m, int * n,              double * A,  int * ldA,  d, e, tauQ,   tauP,             double * work, int * work, int * info )
dgbbrd_ ( vect, int * m, int * n, ncc, kl, ku, Ab, int * ldAb, d, e, Q, int * ldQ, Pt, int * ldPt, C, int * ldC, double * work,        int * info )


QR factorization
geqrf_  ( int * m, int * n, double * A, int * ldA, tau, double * work, int * work, int * info )
Generate Q, after QR
dorgqr_ ( int * m, int * n, k, double * A, int * ldA, tau, double * work, int * work, int * info )
Multiply by Q, after QR, LQ, QL, RQ, TZ factorization
dormqr_ ( side, trans, int * m, int * n, k,    double * A, int * ldA, tau, C, int * ldC, double * work, int * work, int * info )


Reduce to tri-diagonal
dsytrd_ (       uplo, int * n,     double * A, int * ldA,   d, e, tau,    double * work, int * work, int * info )
dsptrd_ (       uplo, int * n,     Ap,       d, e, tau,                 int * info )
dsbtrd_ ( vect, uplo, int * n, kd, Ab, int * ldAb, d, e, Q, int * ldQ, double * work,        int * info )

Generate Q, after tri-diagonalization
dorgtr_ ( uplo, int * n, double * A, int * ldA, tau, double * work, int * work, int * info )
dopgtr_ ( uplo, int * n, Ap, tau, Q, int * ldQ, double * work, int * info )

Multiply by Q, after tri-diagonalization
dormtr_ ( side, uplo, trans, int * m, int * n, double * A, int * ldA, tau, C, int * ldC, double * work, int * work, int * info )
dopmtr_ ( side, uplo, trans, int * m, int * n, Ap,     tau, C, int * ldC, double * work,        int * info )

Eigenvalues of symmetric tri-diagonal
dsterf_ (        int * n, d, e,                                                    int * info )
dsteqr_ ( compz, int * n, d, e, Z, int * ldZ, double * work,                                      int * info )
dpteqr_ ( compz, int * n, d, e, Z, int * ldZ, double * work,                                      int * info )
dstedc_ ( compz, int * n, d, e, Z, int * ldZ, double * work, int * work,                iwork, liwork, int * info )

Eigenvalues of symmetric tri-diagonal (selected)
dstemr_ ( char * jobz, range,        int * n, d, e, vl, vu, il, iu,               int * m,         w, Z, int * ldZ, nzc, isuppz, tryrac, double * work, int * work, iwork, liwork, int * info )
dstegr_ ( char * jobz, range,        int * n, d, e, vl, vu, il, iu, abstol,       int * m,         w, Z, int * ldZ,      isuppz,         double * work, int * work, iwork, liwork, int * info )
dstebz_ (       range, order, int * n,       vl, vu, il, iu, abstol, d, e, int * m, nsplit, w,              iblock, isplit, double * work,        iwork,         int * info )
*/
